{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "id": "DMi9kMlvGFmm",
    "outputId": "4cb57f1e-f403-4668-8b30-5e423790837a"
   },
   "outputs": [],
   "source": [
    "!nvidia-smi"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "ONrKSkSvNquN"
   },
   "source": [
    "# Nova seção"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "5qIf8Ivn6OjQ"
   },
   "source": [
    "# New Section"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "ux-7QMPx6PIV"
   },
   "source": [
    "# New Section"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "n461or5eb8du"
   },
   "source": [
    "# Variational Autoencoders\n",
    "\n",
    "In this notebook, we'll introduce and explore \"variational autoencoders,\" which are a very successful family of models in modern deep learning. In particular we will:\n",
    "\n",
    "\n",
    "1.   Illustrate the connection between autoencoders and classical *Principal Component Analysis (PCA)*\n",
    "3.   Train a non-linear variational auto-encoder that uses a deep neural network\n",
    "\n",
    "### Overview\n",
    "As explained in lecture, variational autoencoders are a way of discovering *latent, low-dimensional structure* in a dataset. In particular, a random data vector $W \\in \\mathbb{R}^d$ can be said to have low-dimensional structure if we can find some functions $e: \\mathbb{R}^d \\to \\mathbb{R}^k$ and $d: \\mathbb{R}^k \\to \\mathbb{R}^d$, with $k \\ll d$, such that $$d(e(W)) \\approx W.$$\n",
    "In other words, $e(W)$ is a parsimonious, $k$-dimensional representation of $W$ that contains all of the information necessary to approximately reconstruct the full vector $W$. Traditionally, $f(W)$ is called an *encoding* of $W$.\n",
    "\n",
    "It turns out that this is meaningless unless we restrict what kinds of functions $d$ and $e$ are allowed to be, because it's possible to write down some (completely ugly) one-to-one function $\\mathbb{R}^d \\to \\mathbb{R}^1$ for any $d$. This gives rise to the notion of *variational autoencoders* where, given some sets of reasonable functions $F$ and $G$, we aim to minimize \n",
    "$$\\mathbb{E}[\\mathrm{loss}(W, d(e(W))]$$ over functions $d \\in D$ and $e \\in E$. As usual, this is done by minimizing the sample analog. \n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "rfey-R_fhGTT"
   },
   "source": [
    "## Linear Autoencoders and PCA: Theory\n",
    "\n",
    "\n",
    "Let $W$ be a centered random vector in $\\mathbb{R}^d$, and let $\\Sigma_n = \\mathbb{E}_n[W_iW_i'] \\in \\mathbb{R}^{d \\times d}$ be its covariance matrix.\n",
    "\n",
    "\n",
    "Consider mutually orthogonal rotations $$X_{ik}:= c_k'W_i,$$ of original $W_i$'s,  where  $$c_\\ell'c_k = 0 \\text{ for } \\ell \\neq k \\text{ and } c_k'c_k=1 \\text{ for each } k.$$ \n",
    "\n",
    "The condition on $c_k's$ ensures that for all $\\ell \\neq k$:\n",
    " $$\n",
    "  \\mathbb{E}_n X_{ik} X_{i\\ell} =0.$$ \n",
    "The rotations $X_{ik}$ are called principal components of $W_i$. \n",
    "\n",
    "In applications, $W_i$ represent high-dimensional raw features (images, for example), and $X^K_{i} = (X_{ik})_{k=1}^K$  represent a lower-dimensional\n",
    "encoding or embedding of $W_i$.   \n",
    "\n",
    "\n",
    "The principal components can be seen as the solution of least squares problem: minimize\n",
    "$$\n",
    "\\sum_{j} \\mathbb{E}_n  (W_{ij} - \\hat W_{ij})^2$$\n",
    "subject to\n",
    "$$\n",
    "\\hat W_{ij} :=  a_j' X^K_{ik}, \\quad X_{ik}:= c_k'W_i.\n",
    "$$\n",
    "Therefore they gives the variational autoencoder corresponding to linear classes with  $\\text{loss}(w, v) = ‖w−v‖^2$ .\n",
    "\n",
    "\n",
    "\n",
    "\\\\\n",
    "\n",
    "\n",
    "**Some Theoretical Details$^*$.** Let $c_1, \\ldots c_d$ denote eigenvectors of $\\Sigma$, where the $c_k$ are normalized so that $\\|c_k\\| = 1$ and listed such that corresponding eigenvalues $\\lambda_k$, which satisfy $\\Sigma c_k = \\lambda_kv_k$, are decreasing in $k$. \n",
    "\n",
    "**Lemma.** The following statements are equivalent.\n",
    "1. Vectors $c_1, c_2, \\ldots c_d$ are the eigenvectors of $\\Sigma$\n",
    "2. For each $1 \\le k \\le d$, the subspace $E_k = \\mathrm{span}(c_1 \\ldots c_k)$ minimizes $$\\mathbb{E}_n\\|W - \\Pi_S W\\|^2$$ over subspaces $S\\subset \\mathbb{R}^d$ of dimension $k$, where $\\Pi_S$ denotes the orthogonal projection onto $S$\n",
    "3. For each $1 \\le k \\le d$, the minimum of $$\\mathbb{E}_n\\|W - ABW\\|^2$$ over matrices $A \\in\\mathbb{R}^{d \\times k}$ and $B \\in \\mathbb{R}^{k \\times d}$ is attained at $A = C^k$, $B=(C^k)'$, where $C^k$ is the $d \\times k$ matrix whose $j^{\\text{th}}$ column is $c_j$. \n",
    "\n",
    "To interpret, since the matrices $A$ and $B$ of statement 3 are the same as linear functions between $\\mathbb{R}^d$ and $\\mathbb{R}^k$, the lemma shows that PCA gives the variational autoencoder corresponding to *linear* classes  with square loss."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "5XzJoBaIvksv"
   },
   "source": [
    "## Linear Autoencoders and PCA: Practice\n",
    "\n",
    "Having just proved that linear autoencoders are the same as PCA, let's do a small sanity check. In particular, let's perform PCA two ways: first using a standard (linear algebra) toolkit, and second as a linear autoencoder using a neural network library. \n",
    "If all goes well, they should give you the same reconstructions! \n",
    "\n",
    "To make it a bit more fun, we will use the [*Labeled Faces in the Wild*](https://www.kaggle.com/jessicali9530/celeba-dataset) dataset which consists of standardized images of roughly 5,000 celebreties' faces. In this data, PCA amounts to looking for a small number of \"proto-faces\" such that a linear combination of them can accurately reconstruct any celebrity's face. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "id": "DQfDiQ8aqpyE",
    "outputId": "d9a71466-ecd7-477f-c508-de08136da7e0"
   },
   "outputs": [],
   "source": [
    "# First, let's download and inspect the data!\n",
    "from sklearn.datasets import fetch_lfw_people\n",
    "faces = fetch_lfw_people(slice_=(slice(70, 198), slice(78, 174)))\n",
    "\n",
    "# 3D Array \"faces.images\" contains images as 2d arrays, stacked along dimension 0\n",
    "n_examples, height, width = faces.images.shape\n",
    "\n",
    "# 2D Array \"design_matrix\" encodes each image as a 1d numeric row, as is conventional in statistics\n",
    "design_matrix = faces.images.reshape((n_examples, -1))\n",
    "\n",
    "n_features = design_matrix.shape[1]\n",
    "\n",
    "print(\n",
    "    \"Labeled Faces in the Wild Dataset: \\n\\\n",
    "    Number of examples: {}\\n\\\n",
    "    Number of features: {}\\n\\\n",
    "    Image height: {}\\n\\\n",
    "    Image width: {}\".format(\n",
    "        n_examples,\n",
    "        n_features,\n",
    "        height,\n",
    "        width))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 326
    },
    "id": "PX_E23v-5yZY",
    "outputId": "2ad8d042-2ea9-4af2-8a19-a71c505c77a1"
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "# Let's gather all the images corresponding to Arnold Scwarzenegger to use as examples\n",
    "\n",
    "# Make a list (of length one!) of labels corresponding to Arnold\n",
    "# Array \"faces.target_names\" tells us which numeric label (index) corresponds to which person name (value)\n",
    "clint_labels = np.where(faces.target_names == 'Clint Eastwood')\n",
    "\n",
    "# Get indices of all images corresponding to this label\n",
    "# Array \"faces.target\" tells us which image (index) corresponds to which numeric image labels (value)\n",
    "clint_pics = np.where(np.isin(faces.target, clint_labels))[0]\n",
    "\n",
    "# Make a helper function so we can take a look at our target images\n",
    "def plot_faces(images, n_row=2, n_col=3):\n",
    "    \"\"\"Helper function to plot a gallery of portraits\"\"\"\n",
    "    plt.figure(figsize=(3.5* n_col, 2.2 * n_row))\n",
    "    plt.subplots_adjust(0.6, 0.5, 1.5, 1.5)\n",
    "    for i in range(n_row * n_col):\n",
    "        plt.subplot(n_row, n_col, i + 1)\n",
    "        plt.imshow(images[i].reshape((height, width)), cmap=plt.cm.gray)\n",
    "        plt.xticks(())\n",
    "        plt.yticks(())\n",
    "    plt.tight_layout()\n",
    "    plt.show()\n",
    "\n",
    "# Let's try it out!\n",
    "plot_faces(\n",
    "    faces.images[clint_pics[:6],:,:] #first six images of Arnold appearing in the dataset\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "Q74vWlBxyIgA"
   },
   "outputs": [],
   "source": [
    "# 1. Find the first 32 principal components of the dataset using the Scikit-learn library \n",
    "# For extra fun, you can do so directly using the singular value decomposition (your mileage may vary!)\n",
    "\n",
    "# We'll use a standard library, which uses linear algebra to compute the principal components.\n",
    "from sklearn.decomposition import PCA\n",
    "\n",
    "# There's no need to de-mean the data. Can you explain why?\n",
    "pca = PCA(n_components=256, svd_solver='randomized').fit(design_matrix)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 326
    },
    "id": "sLi4k8t3DrHe",
    "outputId": "6971e331-fd9d-4f7b-a250-1467f1e6b5e3"
   },
   "outputs": [],
   "source": [
    "# 2. Plot the first 6 \"eigenfaces,\" the six images whose linear span best explains the variation in our dataset\n",
    "eigenfaces = pca.components_\n",
    "plot_faces(eigenfaces[:6])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 326
    },
    "id": "Gmj2lpTfCXKC",
    "outputId": "10366bf3-607a-4ab0-9862-a2ce4d4c86be"
   },
   "outputs": [],
   "source": [
    "# 3. Plot Arnold's face (any image will do!) reconstructed using 1, 8, 64, and 256 principal components\n",
    "face_vector = design_matrix[clint_pics[1]]\n",
    "\n",
    "def reconstruct(image_vector, n_components):\n",
    "  return eigenfaces[:n_components].T @ (eigenfaces[:n_components] @ image_vector.reshape((-1,1)))\n",
    "\n",
    "reconstructions = [reconstruct(face_vector, k) for k in [1,2,4,8,64,256]]\n",
    "plot_faces(reconstructions)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "eoZ_BsXYDE7P"
   },
   "outputs": [],
   "source": [
    "# 4. Train linear autoencoder with 8, 64, and 256 neurons using Keras (example below has dim. 64)\n",
    "# 5. Compare reconstructions of Arnold's face both using MSE and visually"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "urlMaifVJCDc"
   },
   "outputs": [],
   "source": [
    "import tensorflow as tf\n",
    "from tensorflow.keras import Model, Input\n",
    "from tensorflow.keras.layers import Dense\n",
    "\n",
    "encoding_dimension = 256\n",
    "input_image = Input(shape=(n_features,))\n",
    "encoded = Dense(encoding_dimension, activation='linear')(input_image)\n",
    "decoded = Dense(n_features, activation='linear')(encoded)\n",
    "\n",
    "autoencoder = Model(input_image, decoded)\n",
    "\n",
    "autoencoder.compile(optimizer='adam', loss='mse')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "id": "5OTUbWg8NcIE",
    "outputId": "0db2253e-fd71-4875-8ee4-a0f5418a1937"
   },
   "outputs": [],
   "source": [
    "autoencoder.fit(design_matrix, design_matrix,\n",
    "                epochs=100,\n",
    "                batch_size=20,\n",
    "                shuffle=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 185
    },
    "id": "omt7zwGfOOQF",
    "outputId": "14df6ddb-132e-4bb6-aabf-cbfbfc472f97"
   },
   "outputs": [],
   "source": [
    "# Compute neural reconstruction\n",
    "reconstruction = autoencoder.predict(face_vector.reshape(1,-1))\n",
    "\n",
    "# Do visual comparison\n",
    "plot_faces([reconstructions[5],reconstruction],n_row=1,n_col=2)\n",
    "\n",
    "# Do numeric comparison\n",
    "# We also normalize the black/white gradient to take values in [0,1] (divide by 255)\n",
    "print('Mean-squared discrepancy: {}'.format(np.mean(np.power((reconstructions[4].T - reconstruction)/255,2))))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "V6kUoOaWvtxY"
   },
   "source": [
    "Note that NN-programmed auto-encoder is not doing as well as PCAs. Can we improve performance by doing changing parameters (training with more epochs, changing batch size, using penalization)?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "zrdD6A55yJML"
   },
   "source": [
    "## Neural Autoencoders\n",
    "\n",
    "Finally, let's train a nonlinear autoencoder for the same data where $F$ and $G$ are neural networks, and we restrict the dimension to be $k=64$. \n",
    "\n",
    "Visually compare the reconstructions some of Arnold's faces. How much better does the convolutional (nonlinear) model perform than the linear one? \n",
    "\n",
    "The convolutional NNs feature layers that \"summarize\" information by doing local averages over the segments of the immage. (For more details, see, e.g. Smola et al book).\n",
    "\n",
    "**Note: you will want to click ```Edit > Notebook Settings``` and enable a GPU hardware accelerator, otherwise this will take excruciatingly long.**\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "id": "DSN-UliDzXyg",
    "outputId": "8f98a155-018d-42a3-8a1e-f80081c306e3"
   },
   "outputs": [],
   "source": [
    "from tensorflow.keras.layers import Add, Conv2D, Conv2DTranspose, MaxPool2D, Flatten, Reshape, Concatenate, Dropout,GaussianNoise\n",
    "\n",
    "# 1. Train convolutional autoencoder using Keras\n",
    "latent_dim = 256\n",
    "\n",
    "input_image = Input(shape=(height*width,))\n",
    "y = Reshape((height, width, 1))(input_image)\n",
    "y = GaussianNoise(0.1)(y, training=True)\n",
    "x = Conv2D(latent_dim*2, 3, activation=\"relu\",padding='same')(y)\n",
    "x = MaxPool2D(2)(x)\n",
    "x = Conv2D(latent_dim*4, 3, activation=\"relu\",padding='same')(x)\n",
    "x = MaxPool2D(2)(x)\n",
    "x = Conv2D(latent_dim*4, 3, activation=\"relu\",padding='same')(x)\n",
    "x = MaxPool2D(2)(x)\n",
    "x = Flatten()(x)\n",
    "y = Flatten()(y)\n",
    "x = Concatenate()([x,y])\n",
    "encoding = Dense(latent_dim, activation=\"relu\")(x)\n",
    "\n",
    "z = Dense(16*12*latent_dim*4,activation=\"relu\")(encoding)\n",
    "# z = Dense(16*12*latent_dim*4,activation=\"relu\")(z)\n",
    "z = Reshape((16,12,latent_dim*4))(z)\n",
    "z = Conv2DTranspose(latent_dim*2,3,strides=2, activation=\"relu\", padding='same')(z)\n",
    "z = Conv2DTranspose(1,3, strides=2, activation=\"relu\", padding='same')(z)\n",
    "output_image = Reshape((height*width,))(z)\n",
    "convolutional_autoencoder = Model(input_image, output_image)\n",
    "convolutional_autoencoder;\n",
    "convolutional_autoencoder.summary()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "q2TBEiBYeyZD"
   },
   "source": [
    "# 149.5731"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 1000
    },
    "id": "Ulhjit_zxd20",
    "outputId": "c1017900-e634-4e41-8705-e5e320c4f60e"
   },
   "outputs": [],
   "source": [
    "opt = tf.keras.optimizers.Adamax()\n",
    "convolutional_autoencoder.compile(optimizer=\"Adamax\", loss='mse')\n",
    "convolutional_autoencoder.fit(design_matrix, design_matrix,\n",
    "                epochs=300,\n",
    "                batch_size=50,\n",
    "                shuffle=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "X5fxbHKzRLmP"
   },
   "source": [
    "Results\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 185
    },
    "id": "ByufPSRMVSIr",
    "outputId": "3f870955-ddec-453e-cad0-ef8aa7f02454"
   },
   "outputs": [],
   "source": [
    "# Compute neural reconstruction\n",
    "reconstruction = convolutional_autoencoder.predict(face_vector.reshape(1,-1))\n",
    "\n",
    "# Do visual comparison\n",
    "plot_faces([reconstructions[5],reconstruction],n_row=1,n_col=2)\n",
    "print('Mean-squared discrepancy: {}'.format(np.mean(np.power((reconstructions[5].T - reconstruction)/255,2))))"
   ]
  }
 ],
 "metadata": {
  "accelerator": "GPU",
  "colab": {
   "collapsed_sections": [],
   "name": "Autoencoders.ipynb",
   "provenance": []
  },
  "kernelspec": {
   "display_name": "Python 3",
   "name": "python3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}
